home *** CD-ROM | disk | FTP | other *** search
/ Meeting Pearls 1 / Meeting Pearls Vol 1 (1994).iso / installed_progs / gfx / lise2.1 / lise / src / scalc.c < prev    next >
Encoding:
C/C++ Source or Header  |  1993-03-31  |  13.4 KB  |  509 lines
  1. /********************************************************************/
  2. /*                                                                  */
  3. /*    MODUL              : SCALC.C                                  */
  4. /*    VERSION            : 1.0                                      */
  5. /*    DATE               : 24/09/91                                 */
  6. /*    WRITTEN BY         : Christoph Bohle                          */
  7. /*    Modified for spectra calculation by Rainer Kowallik, 06/11/91 */
  8. /* ---------------------------------------------------------------- */
  9. /*    LANGUAGE           : C                                        */
  10. /*    OS                 : UNIX                                     */
  11. /*                                                                  */
  12. /********************************************************************/
  13.  
  14. /********************************************************************/
  15. /*   INCLUDE-Files                                                  */
  16. /********************************************************************/
  17.  
  18. #include <stdio.h>
  19. #include <math.h>
  20. #include <string.h>
  21. #include <ctype.h>
  22. #include <stdlib.h>
  23. #include <spec.h>
  24.  
  25. #ifndef RAND_MAX
  26. #define RAND_MAX 2147483647
  27. #endif
  28.  
  29. int spclenmax = 0;
  30.  
  31. help()
  32. {
  33. printf("scalc  spectra calculator. Call convention:\n");
  34. printf("    scalc [newspc =] formula(function .. number .. spectrum ..)\n");
  35. printf("\nIf no output spectrum is specified, only the first element is\n");
  36. printf("printed on standard output.\n");
  37. printf("Outputspectrum can be specified with a range.\n");
  38. printf("If outputspec allready exists, only the specified part will\n");
  39. printf("be substituted by the result of the calculation.\n");
  40. printf("No calculation will be done with the error and time arrays !\n");
  41. printf("\nPossible function:\n");
  42. printf("exp(expression)\n");
  43. printf("log(expression)\n");
  44. printf("log10(expression)\n");
  45. printf("sin(expression)\n");
  46. printf("cos(expression)\n");
  47. printf("tan(expression)\n");
  48. printf("asin(expression)\n");
  49. printf("acos(expression)\n");
  50. printf("atan(expression)\n");
  51. printf("sinh(expression)\n");
  52. printf("cosh(expression)\n");
  53. printf("tanh(expression)\n");
  54. printf("sqrt(expression)\n");
  55. printf("rnd()                  returns a random number between 0 and 1\n");
  56. printf("lin()                  for(n=1;n<_MAXSPCLEN;n++) spc[n] = n;\n");
  57. printf("sum(spc)               for(n=1;n<_MAXSPCLEN;n++) sum = sum + spc[n];\n");
  58. printf("pi()                   for(n=1;n<_MAXSPCLEN;n++) spc[n] = 3.1415927;\n");
  59. printf("\nOperators are +,-,*,/ and ^\n");
  60. printf("   programmed by Christoph Bohle and Rainer Kowallik\n");
  61. exit(0);
  62. }
  63.  
  64.  
  65.        
  66. /********************************************************************/
  67. /*    Functions                                                     */
  68. /********************************************************************/
  69.  
  70.  
  71. float *expression();
  72.  
  73. char *end_of_brackets( str )
  74. char *str;  
  75. {                                   
  76.   int i;
  77.   i=1;
  78.   do
  79.   {    
  80.     str++;
  81.     switch( *str )
  82.     {
  83.       case '('  : i++;
  84.                   break;
  85.       case ')'  : i--;
  86.                       break;
  87.       case '\0' : i=0;
  88.     }
  89.   }
  90.   while( i>0 );
  91.   return( str );
  92. }
  93.         
  94. float *factor( facstr )
  95. char *facstr;
  96. {
  97. int i,n;
  98. float *result, *tim, *err;
  99. float x,y;
  100. char *help, c;
  101. char comment[80];
  102.   
  103.   /* Float */
  104.   if( isdigit(*facstr) || *facstr=='.' || *facstr=='-' || *facstr=='+' ) {
  105.     result = (float *) calloc(_MAXSPCLEN+2,sizeof(float));
  106.     y = atof(facstr);
  107.     for(i = 0; i < _MAXSPCLEN; i++) result[i] = y;
  108.     return(result);
  109.   }
  110.   
  111.   /* Expression in brackets */
  112.   if( *facstr=='(' )
  113.   {
  114.     *(end_of_brackets( ++facstr )) = '\0';
  115.     return( expression( facstr ) );
  116.   }
  117.   
  118.   /* Function */
  119.   help = facstr;
  120.   while( isalnum(*help) || *help == '_' || *help == '.'
  121.                         || *help == '[' || *help == ':' || *help == ']') help++;
  122.   if( help != facstr )
  123.   {
  124.     c = *help;
  125.     *help = '\0';
  126.     
  127.     /* Function */
  128.     if( strcmp( facstr , "exp" ) == 0 )
  129.     {
  130.       *help = c;
  131.       result = factor( help );
  132.       for(i = 0; i < spclenmax; i++) result[i] = exp(result[i]);
  133.       return( result );
  134.     }
  135.     if( strcmp( facstr , "log" ) == 0 )
  136.     {
  137.       *help = c;
  138.       result = factor( help );
  139.       for(i = 0; i < spclenmax; i++) result[i] = log(result[i]);
  140.       return( result );
  141.     }
  142.     if( strcmp( facstr , "log10" ) == 0 )
  143.     {
  144.       *help = c;
  145.       result = factor( help );
  146.       for(i = 0; i < spclenmax; i++) result[i] = log10(result[i]);
  147.       return( result );
  148.     }
  149.     if( strcmp( facstr , "sin" ) == 0 )
  150.     {
  151.       *help = c;      result = factor( help );
  152.       for(i = 0; i < spclenmax; i++) result[i] = sin(result[i]);
  153.       return( result );
  154.     }
  155.     if( strcmp( facstr , "cos" ) == 0 )
  156.     {
  157.       *help = c;
  158.       result = factor( help );
  159.       for(i = 0; i < spclenmax; i++) result[i] = cos(result[i]);
  160.       return( result );
  161.     }
  162.     if( strcmp( facstr , "tan" ) == 0 )
  163.     {
  164.       *help = c;
  165.       result = factor( help );
  166.       for(i = 0; i < spclenmax; i++) result[i] = tan(result[i]);
  167.       return( result );
  168.     }
  169.     if( strcmp( facstr , "asin" ) == 0 )
  170.     {
  171.       *help = c;
  172.       result = factor( help );
  173.       for(i = 0; i < spclenmax; i++) result[i] = asin(result[i]);
  174.       return( result );
  175.     }
  176.     if( strcmp( facstr , "acos" ) == 0 )
  177.     {
  178.       *help = c;
  179.       result = factor( help );
  180.       for(i = 0; i < spclenmax; i++) result[i] = acos(result[i]);
  181.       return( result );
  182.     }
  183.     if( strcmp( facstr , "atan" ) == 0 )
  184.     {
  185.       *help = c;
  186.       result = factor( help );
  187.       for(i = 0; i < spclenmax; i++) result[i] = atan(result[i]);
  188.       return( result );
  189.     }
  190.     if( strcmp( facstr , "sinh" ) == 0 )
  191.     {
  192.       *help = c;
  193.       result = factor( help );
  194.       for(i = 0; i < spclenmax; i++) result[i] = sinh(result[i]);
  195.       return( result );
  196.     }
  197.     if( strcmp( facstr , "cosh" ) == 0 )
  198.     {
  199.       *help = c;
  200.       result = factor( help );
  201.       for(i = 0; i < spclenmax; i++) result[i] = cosh(result[i]);
  202.       return( result );
  203.     }
  204.     if( strcmp( facstr , "tanh" ) == 0 )
  205.     {
  206.       *help = c;
  207.       result = factor( help );
  208.       for(i = 0; i < spclenmax; i++) result[i] = tanh(result[i]);
  209.       return( result );
  210.     }
  211.     if( strcmp( facstr , "sqrt" ) == 0 )
  212.     {
  213.       *help = c;
  214.       result = factor( help );
  215.       for(i = 0; i < spclenmax; i++) result[i] = sqrt(result[i]);
  216.       return( result );
  217.     }
  218.     if( strcmp( facstr , "rnd" ) == 0 )            /* generate random array */
  219.     {
  220.       *help = c;
  221.       result = (float *) calloc(_MAXSPCLEN+2,sizeof(float));
  222.       for(i = 0; i < _MAXSPCLEN; i++) result[i] = ((float) rand() / (RAND_MAX + 1.0)) ;
  223.       return( result );
  224.     }
  225.     if( strcmp( facstr , "lin" ) == 0 )            /* generate linear array */
  226.     {
  227.       *help = c;
  228.       result = (float *) calloc(_MAXSPCLEN+2,sizeof(float));
  229.       for(i = 0; i < _MAXSPCLEN; i++) result[i] = (float) i ;
  230.       return( result );
  231.     }
  232.     if( strcmp( facstr , "pi" ) == 0 )
  233.     {
  234.       *help = c;
  235.       result = (float *) calloc(_MAXSPCLEN+2,sizeof(float));
  236.       for(i = 0; i < _MAXSPCLEN; i++) result[i] = 3.1415927 ;
  237.       return( result );
  238.     }
  239.     if( strcmp( facstr , "sum" ) == 0 )
  240.     {
  241.       *help = c;
  242.       result = factor( help );
  243.       y = 0.0;
  244.       for(i = 0; i < spclenmax; i++) y = y + result[i];
  245.       for(i = 0; i < spclenmax; i++) result[i] = y;
  246.       return( result );
  247.     }
  248.  
  249. /*  read spectra array */
  250.     result = (float *) calloc(_MAXSPCLEN+2,sizeof(float));
  251.     tim = (float *) calloc(_MAXSPCLEN+2,sizeof(float));
  252.     err = (float *) calloc(_MAXSPCLEN+2,sizeof(float));
  253.     i = readspec(facstr,result,err,tim,comment);
  254.     if(i > spclenmax) spclenmax = i;
  255.     free(err); free(tim);
  256.     return( result );
  257.   }
  258.  
  259.   /* Default */
  260.   if( *facstr == '\0' ) {
  261.     result = (float *) calloc(_MAXSPCLEN+2,sizeof(float));
  262.     for(i = 0; i < _MAXSPCLEN; i++) result[i] = 0.0;
  263.     return( result );
  264.   } else {
  265.     return( factor( facstr+1 ) );
  266.   }
  267. }
  268.  
  269. float *po( postr )
  270. char postr[];
  271. {
  272. char *help;
  273. float *ar1, *ar2;
  274. int i,n;
  275.  
  276.   help = postr;
  277.   while( (*help != '^') && (*help != '\0' ) )
  278.   {
  279.     if( *help == '(' )
  280.         help = end_of_brackets( help );
  281.     help++;
  282.   }
  283.   
  284.   switch( *help )
  285.   {
  286.       case '^' : *help = '\0';
  287.                ar1 = factor( postr );
  288.                ar2 = po( help + 1);
  289.                for(i = 0; i < spclenmax; i++) ar1[i] = pow(ar1[i],ar2[i]);
  290.                free(ar2);
  291.                return( ar1 );
  292.     case '\0' : if( help == postr ) {
  293.                   ar1 = (float *) calloc(_MAXSPCLEN+2,sizeof(float));
  294.                   for(i = 0; i < spclenmax; i++) ar1[i] = 0.0;
  295.                   return( ar1 );
  296.                 } else {
  297.                   return( factor( postr ) );
  298.                 }
  299.   }
  300. }
  301.  
  302. float *term( termstr )
  303. char termstr[];
  304. {
  305. char *help;
  306. float *ar1, *ar2;
  307. int i,n;
  308.   
  309.   help = termstr;
  310.   while( (*help != '*') && (*help != '/') && (*help != '\0' ) )
  311.   {
  312.     if( *help == '(' )
  313.         help = end_of_brackets( help );
  314.     help++;
  315.   }
  316.   
  317.   switch( *help )
  318.   {
  319.      case '*' : *help = '\0';
  320.                ar1 = po( termstr );
  321.                ar2 = term( help + 1);
  322.                for(i = 0; i < spclenmax; i++) ar1[i] = ar1[i] * ar2[i];
  323.                free(ar2);
  324.                return( ar1 );
  325.     case '/' : *help = '\0';
  326.                ar1 = po( termstr );
  327.                ar2 = term( help + 1);
  328.                for(i = 0; i < spclenmax; i++) {
  329.                   if(ar2[i] != 0.0) ar1[i] = ar1[i] / ar2[i];
  330.                }
  331.                free(ar2);
  332.                return( ar1 );
  333.     case '\0' : if( help == termstr ) {
  334.                   ar1 = (float *) calloc(_MAXSPCLEN+2,sizeof(float));
  335.                   for(i = 0; i < spclenmax; i++) ar1[i] = 0.0;
  336.                   return( ar1 );
  337.                 } else {
  338.                   return( po( termstr ) );
  339.                 }
  340.   }
  341. }
  342.  
  343. float *expression( exprstr )
  344. char exprstr[];
  345. {
  346. char *help;
  347. float *ar1, *ar2;
  348. int i,n;
  349.  
  350.   help = exprstr;
  351.   while( (*help != '+') && (*help != '-') && (*help != '\0' ) )
  352.   {
  353.     if( *help == '(' )
  354.         help = end_of_brackets( help );
  355.     help++;
  356.   }
  357.   
  358.   switch( *help )
  359.   {
  360.      case '+' : *help = '\0';
  361.                ar1 = term( exprstr );
  362.                ar2 = expression( help + 1);
  363.                for(i = 0; i < spclenmax; i++) ar1[i] = ar1[i] + ar2[i];
  364.                free(ar2);
  365.                return( ar1 );
  366.     case '-' : *help = '\0';
  367.                ar1 = term( exprstr );
  368.                ar2 = expression( help + 1);
  369.                for(i = 0; i < spclenmax; i++) ar1[i] = ar1[i] - ar2[i];
  370.                free(ar2);
  371.                return( ar1 );
  372.     case '\0' : if( help == exprstr ) {
  373.                   ar1 = (float *) calloc(_MAXSPCLEN+2,sizeof(float));
  374.                   for(i = 0; i < spclenmax; i++) ar1[i] = 0.0;
  375.                   return( ar1 );
  376.                 } else {
  377.                   return( term( exprstr ) );
  378.                 }
  379.   }
  380. }
  381.  
  382. float *analyse_function( funstr )
  383. char *funstr;
  384. {
  385. char *help;
  386. float *result;
  387.   
  388.   help=malloc( strlen(funstr)+3 );
  389.   strcpy( help , funstr );
  390.   result = expression( help );
  391.   free( help );
  392.   return( result );  
  393. }
  394.  
  395. /********************************************************************/
  396. /* Main-function                                                    */
  397. /********************************************************************/
  398.  
  399. main( argc , argv  )
  400. int argc;
  401. char *argv[];              
  402. {
  403. int i, l, ix1, ix2;
  404. char *fun, *snam;
  405. float *result, *oldspc, *err;
  406.   
  407.   snam = (char *) malloc(80);
  408.  
  409.   i = checkopt(argc,argv,"-sowatjibtsnich",snam); /* we must call checkopt once */
  410.  
  411.   l = 1;
  412.   for( i=1 ; i < argc ; i++ )
  413.   {
  414.     l += strlen( argv[i] ) + 1;
  415.   }
  416.   fun = malloc( l + 3 ); *fun = '\0';
  417.   for( i=1 ; i < argc ; i++ )
  418.   {
  419.     strcat( fun , argv[i] );
  420.     strcat( fun , " ");
  421.   }
  422.  
  423.   ix1 = 0;
  424.   ix2 = 0;
  425.   oldspc = (float *) calloc(_MAXSPCLEN+2,sizeof(float));
  426.   err = (float *) calloc(_MAXSPCLEN+2,sizeof(float));
  427.   find_res_spc(fun,snam,&ix1,&ix2,oldspc,err);
  428.   result = analyse_function( fun );
  429.  
  430.   if((ix1 != 0) || (ix2 != 0)) {
  431.      if(ix2 < 0) ix2 = spclenmax;
  432.      for(i = ix1; i <= ix2; i++) oldspc[i] = result[i];
  433.      writespec(snam,oldspc,err,spclenmax,2,"test");
  434.      free(err); free(oldspc);
  435.   } else {
  436.      printf("%f\n",(double)result[0]);
  437.   }
  438.   free(snam); free(result);
  439.   exit(0);
  440. }
  441.  
  442. find_res_spc(str,snam,ix1,ix2,oldspc,err)
  443. char *str, *snam;
  444. int  *ix1, *ix2;
  445. float *oldspc, *err;
  446. {
  447. int i,j,l,n;
  448. char c;
  449. char *z, *s;
  450. FILE *stream;
  451. float *tim;
  452.  
  453.  
  454.    z = (char *) malloc(80);
  455.    s = (char *) malloc(80);
  456.  
  457.    l = instr("=",str);
  458.    if(l < 0) return;
  459.  
  460.    snam[0] = 0; i = 1; j = 0; c = str[0];
  461.    while((c != 0) && (c != '[') && (c != '=')) {
  462.       if(c != ' ') snam[j++] = c;
  463.       c = str[i++];
  464.    }
  465.    snam[j] = 0;
  466.  
  467.    if(c == '[') {
  468.       s[0] = 0; j = 0; c = str[i++];
  469.       while((c != 0) && (c != ':') && (c != ']')) {
  470.          if(c != ' ') s[j++] = c;
  471.          c = str[i++];
  472.       }
  473.       s[j] = 0; *ix1 = atoi(s); *ix2 = *ix1;
  474.       if(c == ':') {
  475.          s[0] = 0; j = 0; c = str[i++];
  476.          while((c != 0) && (c != ']')) {
  477.             if(c != ' ') s[j++] = c;
  478.             c = str[i++];
  479.          }
  480.          s[j] = 0; *ix2 = atoi(s);
  481.       }
  482.    } else {
  483.       *ix1 = 0;
  484.       *ix2 = -1;
  485.    }
  486.  
  487.    j = 0; n = l; l = strlen(str);
  488.    for(i = n; i <= l; i++) str[j++] = str[i];
  489.   
  490.    strcpy(z,snam);
  491.    stream = fopen(z,"r");
  492.    if(stream == NULL) {
  493.       strcat(z,".spc");
  494.       stream=fopen(z,"r");
  495.    }
  496.    if(stream!=NULL) {
  497.       fclose(stream);
  498.       tim = (float *) calloc(_MAXSPCLEN+2,sizeof(float));
  499.       i = readspec(snam,oldspc,err,tim,s);
  500.       free(tim);
  501.       if(i > spclenmax) spclenmax = i;
  502.    }
  503.    
  504.    if(*ix2 > spclenmax) spclenmax = *ix2;
  505.    free(s); free(z);
  506.    return(0);
  507. }
  508.  
  509.